squirrels %>% 
  ggplot()+
  geom_bin2d(aes(y=lat,x=long))+
  theme_minimal()+
  coord_map()+
  labs(
    x=NULL,y=NULL,
    fill="Count of\nSquirrels",subtitle="Heat Map of Squirrel Sightings in Central Park",title='Squirrel Town')+
  theme(
    legend.position = 'right',
    plot.caption = element_text(hjust = 0)
    )

plot_data<-final_dataset %>% 
  filter(circuit_label_type_desc=="Stove/Oven/Range") %>% 
  filter(doy<350) %>% 
  group_by(date,date_2020) %>% 
  summarise(kwh=sum(kwh)) %>% 
  arrange(date) %>% 
  group_by(year(date)) %>% 
  mutate(sma=zoo::rollmean(kwh,30,na.pad = T,align = 'right'))

p2<-ggplot(plot_data)+
  geom_smooth(aes(x=date_2020,y=kwh,color=factor(year(date))),linewidth=1.5,method = 'loess',span=.17,se=F)+
  geom_hline(yintercept = 0,color=NA)+
  geom_vline(xintercept = as.Date("2020-03-08"))+
  labs(x="Day of the Year",y="Smoothed Daily Energy Usage (kWh)",color="Time Period",
       title="Comparison of Stove Energy Usage by Day of Year", subtitle = "Stove usage increased dramatically when the Pandemic started")+
  scale_color_brewer(type = "qual",palette =  6)+
  scale_x_date(date_breaks = "1 month",date_labels = "%b",limits = c(as.Date("2020-01-01"),as.Date("2020-09-01")))+
  theme_minimal()+
  theme(
    legend.position = 'bottom',
    panel.grid.minor.x = element_blank()
    )
p2

agg_dat<-wwc_outcomes %>% 
  group_by(team,round) %>% 
  summarise(n=n(),mean=mean(score)) %>% 
  filter(round!="Third Place Playoff") %>% 
  filter("Final"%in%round) %>% 
  group_by(team) %>% 
  filter(n[round=="Final"]>1)

agg_dat$round<-factor(agg_dat$round,levels = c("Group","Round of 16","Quarter Final","Semi Final","Final"))

h2h<-wwc_outcomes %>% 
  filter(team%in%agg_dat$team) %>% 
  group_by(year,yearly_game_id) %>% 
  filter(n()==2)

output<-NULL
for(i in unique(h2h$team)){
  this_summary<-h2h %>% 
    group_by(year,yearly_game_id) %>% 
    filter(i%in%team) %>%
    group_by(this_team=i,team) %>% 
    summarise(
      w=sum(win_status=="Won"),
      l=sum(win_status=="Lost"),
      d=sum(win_status=="Tie")) %>% 
    mutate(
      pct=l/(w+l)
    )
  output<-bind_rows(
    this_summary,
    output
  )
  
}

final_summary<-wwc_outcomes %>% 
  filter(round=="Final") %>% 
  group_by(team) %>% 
  summarise(apps=n(),first=min(year),goals=sum(score)) %>% 
  arrange(-goals)

p3<-output %>% 
  left_join(codes) %>% 
  filter(this_team!=team) %>% 
  ggplot()+
  geom_col(aes(x=country,y=pct,fill=country))+
  theme_minimal()+
  theme(
    legend.position = 'bottom',
    strip.background = element_rect(fill="gray80"),
    strip.placement = 'left',
    strip.text.y.left = element_text(angle = 0))+
  scale_y_continuous(labels = scales::percent,expand = c(0,0,.05,0))+
  scale_fill_manual(values=rev(c("black","#00205B","#BC002D","#FFCC00")))+
  coord_flip()+
  facet_grid(paste(this_team,"V.")~.,scale="free_y",switch = "y")+
  guides(fill='none')+
  labs(x=NULL,y="Winning Percentage (excluding draws)",
       title="Better v. Best",subtitle = "The United States Women's National Soccer team wins 75% of its matches against other top teams")

p3

sites_with_mains<-subset(points,circuit_label_type_desc=="Mains")$ee_site_id
these_eus<-c(
  "Mains","Hot Tub","Electric Resistance Storage Water Heaters",
  "Central AC", "Clothes Washer", "Clothes Dryer","Refrigerator/Freezer",
  "Other","Stove/Oven/Range","Electric Vehicle Charger","Dishwasher",
  "Ductless Heatpump","Electric Furnace","Ducted Heatpump"
  )

#Step 1 drop incomplete circuits
circuit_agg<-analysis_data %>% 
  filter(ee_site_id%in%sites_with_mains) %>% 
  group_by(ee_site_id,regname) %>% 
  summarise(n=n())

retain_circuits<-circuit_agg %>% 
  filter(n>12755)

analysis_data<-analysis_data %>% 
  inner_join(retain_circuits %>% select(-n))

#Step 3 drop errant days
day_agg<-analysis_data %>% 
  group_by(date,hour) %>% 
  summarise(n=n())

retain_hours<-day_agg %>% 
  filter(n>625) %>% 
  filter(!date%in%as.Date(c("2019-09-25","2019-09-26")))

analysis_data<-analysis_data %>% 
  inner_join(retain_hours %>% select(-n))

weather_stations <- read_csv("data/NEEA HEMS/weather_stations_geocodio.csv",name_repair = 'universal') %>% 
  mutate(station_id=paste(USAFID,WBAN,sep="-")) %>% 
  mutate(State=case_when(
    State=="OR"~"Oregon",
    State=="WA"~"Washington",
    State=="ID"~"Idaho",
    State=="MT"~"Montana"
  ),
  County=gsub(" County","",County,fixed = T))

state_maps<-map_data("state",c("oregon","washington","idaho","montana"))
county_maps<-map_data("county",c("oregon","washington","idaho","montana"))

sites_ws<-sites %>% 
  filter(ee_site_id%in%analysis_data$ee_site_id) %>% 
  left_join(
    weather_stations %>% 
      group_by(station_id) %>% 
      filter(1:n()==1)
  ) 

data_for_map<-analysis_data %>% 
  left_join(sites_ws %>% select(ee_site_id,County,State,LAT,LONG)) %>% 
  group_by(County,State) %>% 
  summarise(sites=n_distinct(ee_site_id), data=n(),LAT=first(LAT),LONG=first(LONG))

county_maps<-county_maps %>% 
  mutate(across(region:subregion,stringr::str_to_title)) %>% 
  left_join(data_for_map,by=c("region"="State","subregion"="County"))

city_labs<-maps::us.cities %>% 
  filter(country.etc%in%c("OR","ID","WA","MT")) %>% 
  filter(name%in%c("Seattle WA","Bend OR","Portland OR","Boise ID","Spokane WA","Missoula MT")) %>% 
  mutate(lat=ifelse(name%in%c("Bend OR"),lat-.5,lat+.5),
         long=ifelse(name%in%c("Seattle WA"),long-2,long))

p1<-ggplot()+
  geom_polygon(data=state_maps,aes(x=long,y=lat,group=group),color='black',fill=NA)+
  geom_polygon_interactive(data=county_maps,aes(x=long,y=lat,group=group,fill=data,tooltip = paste(subregion,region), data_id = paste(subregion,region)),color=NA,size=.1)+
  # geom_point(data=data_for_map,aes(x=LONG,y=LAT,size=data,fill=sites),color="black",shape=21,alpha=.5)+
  geom_label(data=city_labs,aes(x=long,y=lat,label=name))+
  coord_map(xlim = c(-127,-110))+
  theme_void()+
  scale_size_continuous(labels = scales::comma)+
  scale_fill_viridis_c(na.value = "white",labels=scales::comma)+
  theme(
    legend.key.width=unit(.4,"in"),
    legend.position = 'bottom',
    legend.box= 'horizontal')+
  labs(fill="Number of Observations",size="Total Observations",
       title="Geographic Distrobution of NEEA HEMS Data")+
   guides(fill = guide_colourbar(title.position="top", title.hjust = 0.5),
         size = guide_legend(title.position="top", title.hjust = 0.5))

eu_summary<-analysis_data %>% 
  left_join(points %>% select(ee_site_id,regname,circuit_label_type_desc)) %>% 
  left_join(sites_ws %>% select(ee_site_id,region=State,subregion=County)) %>% 
  group_by(region,subregion,circuit_label_type_desc) %>% 
  summarise(
    n=n()
  ) %>% 
  filter(!circuit_label_type_desc%in%c("Other","Mains","Gas Furnace (Component)")) %>% 
  mutate(
    category=case_when(
      circuit_label_type_desc%in%c("Central AC","Ducted Heatpump","Ductless Heatpump","Electric Baseboard Heaters","Other Zonal Heat","Electric Furnace") ~ "HVAC",
      grepl("water heater",circuit_label_type_desc,ignore.case = T) ~ "Water\nHeating",
      circuit_label_type_desc%in%c("Stove/Oven/Range","Microwave","Garbage Disposal","Dishwasher","Refrigerator/Freezer") ~ "Kitchen",
      TRUE~"Other"
    )
  )

p2<-ggplot(eu_summary)+
  geom_col_interactive(aes(x=n,y=forcats::fct_reorder(abbreviate(circuit_label_type_desc,30),n),fill=category,tooltip = paste(subregion,region), data_id = paste(subregion,region)))+
  scale_x_continuous(labels = scales::comma,expand=c(0,0,.05,0))+
  scale_fill_brewer(type="qual",palette = 2)+
  theme_classic()+
  facet_grid(category~.,scale="free_y",space="free")+
  labs(x="Number of Observations",y="Equipment Type",fill="Equipment Category",
       title="Appliance Distrobution of NEEA HEMS Data",
       subtitle = "Data available for range of household appliances")+
  theme(
    legend.position = 'bottom',
    strip.text.y = element_text(angle = 0),
    title = element_text(size=9))+
  guides(fill='none')

a <- girafe(code = print(p1 + p2 + plot_layout(widths = c(.65,.35))))

a